Finding apparent horizons in numerical relativity.

We review various algorithms for finding apparent horizons in 3+1 numerical relativity. We then focus on one particular algorithm, in which we pose the apparent horizon equation H ≡ ∇in + Kijnn − K = 0 as a nonlinear elliptic (boundary-value) PDE on angular-coordinate space for the horizon shape function r = h(θ, φ), finite difference this PDE, and use Newton’s method or a variant to solve the finite difference equations. We describe a method for computing the Jacobian matrix of the finite differenced H(h) function H(h) by symbolically differentiating the finite difference equations, giving the Jacobian elements directly in terms of the finite difference molecule coefficients used in computing H(h). Assuming the finite differencing scheme commutes with linearization, we show how the Jacobian elements may be computed by first linearizing the continuum H(h) equations, then finite differencing the linearized continuum equations. (This is essentially just the “Jacobian part” of the Newton-Kantorovich method for solving nonlinear PDEs.) We tabulate the resulting Jacobian coefficients for a number of different H(h) and Jacobian computation schemes. We find this symbolic differentiation method of computing the H(h) Jacobian to bemuch more efficient than the usual numerical-perturbation method, and also much easier to implement than is commonly thought. When solving the discrete H(h) = 0 equations, we find that Newton’s method generally shows robust convergence. However, we find that it has a small (poor) radius of convergence if the initial guess for the horizon position contains significant high-spatial-frequency error components, i.e. angular Fourier components varying as (say) cosmθ with m ∼ 8. (Such components occur naturally if spacetime contains significant amounts of high-frequency gravitational radiation.) We show that this poor convergence behavior is not an artifact of insufficient resolution in the finite difference grid; rather, it